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^ ■ Abstract 

. We present a mesoscopic model, based on the Boltzmann Equation, for the 

\^ I interaction between a solid wall and a non-ideal fluid. We present an analytic 

' derivation of the contact angle in terms of the surface tension between the 

liquid-gas, the liquid-solid and the gas-solid phases. We study the dependency 
of the contact angle on the two free parameters of the model, which determine 
the interaction between the fluid and the boundaries, i.e. the equivalent of the 
■ wall density and of the wall-fluid potential in Molecular Dynamics studies. 

^ . We compare the analytical results obtained in the hydrodynamical limit 

| for the density profile and for the surface tension expression with the numer- 

ical simulations. W e compare also our two -phas e appr o ach wi th some exact 



results obtained by lLauga &, Stone! (200I) and IPhilirJ (|l972al lbT) for a pure 



hydrodynamical incompressible fluid based on Navier-Stokes equations with 
boundary conditions made up of alternating slip and no-slip strips. Finally, 
we show how to overcome some theoretic al limitations connec ted with the dis- 
cretized Boltzmann scheme proposed by Shan k, Chenl ( 19941 ) and we discuss 
the equivalence between the surface tension defined in terms of the mechanical 
equilibrium and in terms of the Maxwell construction. 



1 Introduction 



The physics of molecular interactions at fluid-solid interfaces is a very active research 
area with a significant impact on many emerging applicatio ns in material science, 
chemi st ry, micro /nanoeng i neering, biology a nd medicine, see IWhitesides fc Stroock 
fenmh : iGad-eFHakl (|l99QMrTo fc Tail (|l998h . Many problems require the spreading 
of a liquid on a solid that may either be a simple flat and clean surface or present 
some degree of roughness contaminated bv compounds with different chemical- 



effects, ffuid-solid interactions become particularly important for micro- and nano- 
devices, whose p h ysical behavior is largely affec ted by high surface/ volume ratios 
( Tabelind (12003! ) : Karniadakis fc Beskok I l|2002h ) whose direct conseq uence is the 
enhan c ement of capillary ph e nome na with respect to bulk properties (jDe Gennesl 
iRowlinson fc Widoml (|l982| )). On the theoretical side, very little is known 



because of the difficulties to match the classical infinite-volume thermodynamics 
description with surface effects. From the experimental side the study of the surface 



case (IBico et al. 


1999^: Craig et al. 


(12001): Ou et al. (2004): 


Onda et al. 


n 


996 


): Pit et oil (12000 




Vinoexadova & Yabukovl 


(12001. 2002 


); 


Cheng & Giordano! ( 


20021): IChoi et al. (20031)). 



Maurer et al\ (120041): 



20031 ) : IZhufc Granick 



In most cases, to reach quantitative results on specific problems, one is forced 
to rely on numerical simulations, especially in presence of complex boundary condi- 
tions. To date, two major approaches dominate this field from the numerical side. 
The first one is based on a pure hydrodynamical description, with the interaction 
between the flow and the solid fully r enormalized in terms of ad h o c boundary con di- 
tions for the hydrodynamical fields dCottin-Bizonne et al. ( 2004 ) ; Philip! ( 1972a f]) ] 
Prieziev et al\ (12005! ): lLauga fc Stone! (j2003| )). The main drawback is represented 
by the difficulty to describe a variety of different solid properties, with complex 
roughness landscape and chemical-physical attributes. The main advantage is that 
one can directly focus on spatial and frequency variations up to the typical hydro- 
dynamical scales. 

The second approach attacks the problem from an atomistic description, by 
integrating the Newton equations for a set of molecules interacting via a Lennard- 
Jones potential. This is the basic idea behind Molecular Dynamics (MD), which 
requires an additional ad hoc tuning of the free parameters enteri ng the potential be- 
tween li quid-liquid mo l ecules a nd between liquid-sol id molecules dBocouet fc Barratl 
fll99r " 
1990? 



juid-liquid mo l ecules a nd between liquid-sol id molecules llrsocouet .narratl 
Cieolak et all (l200l! ) : IPrieziev et all ()2005!) : iThompson fc Robbins! fll989l 
Thompson fc Troian ( 1997f )). These parameters are fine tuned by compari- 



son with the experiments and are mainly of two types: the overall strength of the 
interactions and the typical interacting distance (fixed by the relative weight be- 
tween the attractive and repulsive terms). The main drawback is here represented 
by the c ongenital scale separat i on be tween this method and continuum phenom- 



Brenner and Ganesan (j2000| )) and consequently the inability to describe 



ena (see 

spatial fluctuations on scales which are much larger than the inter-molecular inter- 
actions and temporal fluctuations larger than a few milliseconds. The advantage 
is given by its apparent " ab-initio" nature, although to be of any practical use, 
the m e thod needs to be supplemented with experimental data (jKoplik fe Banavar 
ril99l! ): lp3^aportl (|l995l) l. 

In this paper we follow a third route, focusing on a mesosco pic modeling of th e 
solid-liquid interaction based on Kinetic Theory of dense fluids ( Cercignanj ( 199lh : 



Chapman and Cowling! ([197 1 



The method, known as Lattice Boltzmann equation (LBE), directly accesses spa- 
tial and temporal fluctuations at the hydrodynamical level with the extra bonus 
of a large flexibility in the des cription of the chemical and physical properties 
of the boundary conditions (see McNamara fc Zanetti ( 1998!): Higuera fc Jimenez! 



Benzi et all (|l992h ). With respect to the MD approach, one pays the price to move 
the matching with the experimental data to the level of bulk interaction between 
Boltzmann distribution function (the equivalent of liquid-liquid MD potential) and 
to the boundary conditions imposed on the Boltzmann equation (the equivalent of 
the solid-li quid potential i n the MP). Ideall y one would supply the infinite BBGKY 
hierarchy (jMohline (jl982h : lKirkwood(jl95(t ) typical of any kinetic description, with 
atomistic information, thereby closing the problem without any approximation. In 
most cases, a more practical approach is taken: in order to derive useful kinetic de- 
scription, some educated guess on the many-body BBGKY hierarchy are proposed 
and tested a-posteriori. 

In this paper, we shall focus on surface effects in presence of phase coexistence 
between a liquid and its saturated vapor. In particular, we aim at investigating how 
to develop an effective mesoscopic description of the surface tension between the 
liquid-gas, ji g , the li quid-solid, 7z.s the gas- s olid -x ^ phases and , more practically of 
the contact ang le 9 kowlinson fc Widoml (|l982f ): IPe Gennesl lj200ah ) that can be 
defined from the above surface tensions: 



cos (9) 



Igs lis 

Tig 



(1) 



We will perf orm an analyti c al and num erical study within the mean field method 
proposed by Shan fc Chen ( 19931 1994 ). based on a Lattice Boltzmann equation 
with an effective two-body potential described only in terms of local, single molecule, 
properties of the fluid (see next section for a detailed description of the method). The 
model provides, to our opinion, the simplest coherent description of the many-body 
interaction typical of dense fluids wi t hin the Latt i ce Boltzmann Equa ti on framework 
for no n -ideal fluids (|Shan fc Chenl (|l998l Il994h : ISwift et all f)l995h : iHe fc Poolen 



f)2nn2h : lKwol3 l|2004 ) 



In this paper, we review first the method as defined bv lShan fc Chel f)l99.4ll99l 
for bulk flows (no boundaries) and we extend it to include the interaction with a 
given solid surface by the introduction of suitable boundary conditions. This defines 
a theoretical scheme able to incorporate non-ideal effects (phase transitions) trig- 
gered by the presence of the solid boundary and complex fluid properties connected 
to the actual density profiles (contact lines, contact angle, capillary phenomena, 
surface tensions, etc.). The main result of the first part is an exact analytical 
expression of the contact angle in terms of the surface tensions derived from the 
hydrodynamical limit of the Boltzmann equations. In the second part, we perform 
a systematic study of the contact angle dependency on the boundary properties 
and we compare our two-phase approach with some exact results obtained for a 
pure hydrodynamical singl e-phase fluid based on Navier-Stokes e quations with suit- 
able boundary conditions ([Lauga fc Stonel (|2003h : IPhilfnl (|l972al lbh. We also show 
how to overcome some theoretical limitat ions connected with the discretized Boltz- 
mann scheme here utilized, proposed in Shan fc Chen ( 1994 ). and we discuss the 
equivalence between the surface tension defined in terms of the me chanical e quilib - 
rium o r in ter ms of the thermodynamical "Maxwell construction" ( Stanley ( 197l| ): 
Huang ( 1987t )). Finally, we discuss the possible application of this method to de- 
scribe non-stationary flows in micro-channel, the apparent slip phenomenon and the 



2 The Shan-Chen (SC) approach to non-ideal flu- 
ids: the inclusion of boundaries. 



As soon as one goes beyond the "ideal gas" description, allowing for phase tran 
sition and for density and temperature variations inside the flow, the c orrect way 



to app roach the kinetic problem is to start from the BBGKY formalism (jKirkwood 



(195O0). Phase transitions are triggered by critical dependency of the thermody- 



namic variables on small variations in the local density, temperature and pressure 
fields. In order to describe such phenomena, one needs to go beyond the description 
based on the probability density to observe a single molecule with a given velocity, 
Vi, at position 7*1 and at time t, fi(ri,Vi,t). In particular, one needs to consider 
at least the two-particle distribution, f\ t2 (ri, Vi, r 2 , v 2 , t), which explicitly enters in 
the Boltzmann equation via the collisional term, Q: 

dth + v 1 -d n h + K x -d vl f 1 = Q, (2) 
where Ki is an external body force and 

tt = - J dv 2 dr 2 d Vl f 1>2 d ri V(r 12 ) (3) 

with V(r 12 ) = V^dri — r 2 |) being the inter-particle potential. The BBGKY hierarchy 
prescribes the evolution of / lj2 in terms of the three particle densities, /i j2j 3, the 
evolution of /i j2| 3 in terms of the four-particles density and so on. The simplest 
closure which takes into account the two-particles interaction consists in adopting a 
"mean field" approach for the collisional terms. This approach starts by rewriting 
/i >2 in the equivalent form: 

/i,2(ri, vi, r 2 , v 2 ) = /i(ri, vi)f 2 (r 2 , v 2 )g(r 1 , r 2 , v 1 , v 2 ) 

where we have introduced the two particles correlation function, g(ri, r 2 ,v±, v 2 ). To 
proceed further one needs to make some approximation on the two-body correlation 
function g(r±, r 2 , vi, v 2 ). The celebrated "molecular chaos" assumptio n of Boltz 



mann gives g = 1, i.e. absence of both velocity and spatial correlations ([Cercignani 



f 199l| )). In the less restrictive case where only velocity correlations vanish, one has: 



g(r 1 ,r 2 ,v 1 ,v 2 ) = g{r x ,r 2 ) 
and the collisional term can be (see also iMartvs! \ 19991 )) rewritten as: 

« = -d v Ji J dr 2 p(r 2 )g(r h r 2 )d ri V(r 12 ) (4) 

where we have used the definition of the local density as 

p(r,t) = J dvf(r,v,t). (5) 



The approximation (HI) is at the core of manv Lattice Boltzmann description of 



Q oc K\-d vl fi, and may be seen as a renormalization of the local pressure tensor vi a 
the intro d uctio n of non-ideal terms in the equation of state ( Chen fc Doolen ( 199§[ ): 
He et gfj(|2002h ). 

More quantitatively, if we start from equation (j2J) and together with (j3J) we define 
the local momentum as 



pu(r,t) = J dvf(r,v,t)v (6) 

we obtain (see iMartvsl dl99i : iKirkwoodl ljl95(t ) conservative equations for the two 
local fields: 

d t p + V • (pu) = (7) 



d t (pu) + V • (J dvfvv) = J dvVtv = -VW 



(8) 



where the tensor W is directly related to the interaction potential V(rv2)'- 

Wij = ~ J dr l dvidv 2 dr 2 fif2g5(r- n) V"(ri2)rf 2 1 (ri 2 )i(ri 2 )j. (9) 

The previous equations a re clearly locally c onservative for density and globally con- 
servative for momentum ( Kirkwoodl ( 195(1 ) ). 

In the realm of Boltzmann Equations, a popular way to simplify the collisional 
integral is to write it as a simple relaxation term (with characteristic time r) towards 
a suitable local equilib rium. This is the celebrated BGK approximation given by 
Bhatnagar et all (jl954j ) and one may wonder which is the simplest BGK description 
consistent with equations (|7l8j) . By consistent BGK description, we mean a single 
time relaxation collisional term of the form: 



n 



BGK 



(/ " / (M) ) 



(10) 



where the equilibrium distribution /( M ) is the local maxwellian in D dimensions: 

.^2 



/ 



(M) 



P 



{2tt KT') D I' 



cxp 



(V — u 
2KT 



and the parameters p'(r,t), u'(r,t), T'(r,t) can be chosen in such a way to be 
consistent with global balance equations. In particular straightforward Gaussian 
integration yields to p'(r,t) = p(r,t) and 



u'(r,t) = u(r,t) - —VW(r,t) 
rp 



(12) 



which means a space-ti me dependent shift of the me an- value of the local Maxwellian. 
The Shan-Chen model dShan fc Chenl (|l99.1 Il994 )) is precisely equivalent to such 
an approach with the assumption 



T7 T/T/ 



( m +\ _ r>. I Ar, \ 



( 1 Q^ 



where with <El v we denote the integration over the angular dependency of the unit 
vector, v and where ip(r,t) is a mean field potential which depends on r and t only 
through the local density, i.e. ip(r,t) = ip(p(r,t)). In the above, (sma^l represents 
the range of the interactions and Qj, a coupling constant, something like an inverse 
temperature for the model. 

In principle, the temperature T'(r, t) in should also be changed on account of 
thermodynamics consistency and total energy conserving dynamics (potential plus 
kinetic). However, being interested in isothermal phenomena, we keep KT' = c^, (? s 
being the sound speed veloci ty (for a possib le extension of the model to include also 
temperature fluctuations see 

Upon discretization of the approximati on described before, we derive imm edi- 
ately the (Lattice) Boltzmann Equations (|Succi ( 2001 ); Wolf-Gladrow ( 2000h ) as 
follows: 



f a (x + dAt, t + At)- f a (x, t) 



■— [fa(x,t)-f ( 
T L 



a 



(p(x,t),u'(x,t))] (14) 



where x runs on a two (or three) dimensional Lattice and At = 1 is the time 
stepping in the numerical scheme. The LHS of (|14j) is the molecular free streeming 
of a discrete set (c a , a = 0, N) of velocities whereas the right-hand side represents 
molecular collisions via a simple relaxation towards the local equilibrium fa 9 ^ (the 
local Maxwellian expanded to second order in the Mach number) in a time lapse 
of the order of r. This relaxation time fixes the fluid kinematic vis c osity as v 



(r - 1/2) where c s = l/VZ in the present work fsee IWolf-Gladrowl (j2000h . The 



fluid density and momentum are given by 



5> 



u 



and can be shown to evolve according to the Navier-Stokes equations of fluid- 
dynamics, Wolf-Gladrow ( 2000j ): 



p[d t u + (u ■ V)u] = 
d t p+V ■ (pu) = 



-VP + F + V- {ypVu) 



(15) 



being P the ideal pressure tensor given by the perfect-gas equation of state: Po = 
I cj,p, where / is the unit tensor. Non-ideal effects are modeled through the 
self-consistent body force term (|13|) that can be discretized as: 



F(x, t) = —G b ip(x, t) 22 w a i){x + c a At, t)c 



(16) 



where ip(x,t) is the lattice version of the mean field potential previously used and 
w a are normalization weights (see Appendix A for more technical details). Due to 
this body force at each time step we consistently re-define the velocity u' in the 
equilibrium distribution as: 



which in turns implies a non-diagonal pressure tensor P deduced from the condition: 



-VP - VP = F. (17) 

By expanding in Taylor series the inter-particle potential one gets (see Appendix A) 
in the case At = 1: 



Si 



'-ctGbdii/jdjip. 



Let us notice that there is always a certain degree of arbitra riness in deriving the 
full ex pression of the pressure tensor in the continuum limit Rowlinson fc WidomI 
( 1982I ). The most reasonable way to do it is to impose that the constraint (|17|) 
is verified up to the second order in the Taylor expansion of both pressure and 
forcing expression. We stop at second order because Navier-Stokes equations are 
obtained from the LBE at the second order in the Chapman Enskog expansion. 
The above expre ssion for the pressur e tensor is different fro m the one prop o sed in 
equation (19) of IShan fc Chenl ( 1994 ). The reason is that in Shan fc Chen ( 1994 ) 
the constraint (J 17)) is verified only up to the first order in Taylor expansion. This 
difference is important for the thermod ynamic consistency o f isothermal flow as it 
will be described in section |2~H see also lHe fc Pooler] ((2002) ■ 

The first two terms in the diagonal part of ()18j) describe the bulk homogeneous 
phase transition by the non-ideal equation of state: 



-p/po ^ 



(20) 



with a typically used functional form: 

Hp) = (i - e 

with p a reference density. 

Now, we will consider a general background that is independent of the functional 
forms of ip, given for granted that, for a quantitative agreement with MD and 
experiments, one can choose different functional forms and/or different inter-particle 
interactions. We will be back on the importance of the functional form later on in 
section ()2.4|) . where we investigate the Maxwell construction in the (P, V) diagram 
of the model. Let us only discuss for the moment the qualitative behavior imposed 
on the p, dependency of ijj(p) by two physical constraints. First, for small densities 
p we need to recover the equation of state of an ideal gas, which requires that: 

tfj(p) oc p p^O. 

Second, for large local densities, the interacting potential must saturate: 

V>(p) — > const. p > p 

a requirement meant to mimic the hard-core properties of real molecules which pre- 
vent unphysical density accumulations. Any smooth functional form which satisfies 
the two previous requirements leads to a phase transition as soon as Gb becomes 
smaller than a critical value Q c (with Q c being negative), where the fluid start to ex- 
hibit two coexisting phases with the same pressure. All the other terms in eq. ()18|) . 
which depend explicitly on the density variations, describe the development of an 
interface profile with its own surface tension (as soon as the isotropy of the pressure 



2.1 SC model without boundaries 



For the sake of completeness, let us summarize again the steps needed to calculate 
the density profile in prese nce of a liquid - gas in terface in an unbounded domain as 
shown for the first time bv lShan Sz Chenl (J1994J). This calculation will be used later 
on to implement the expression of the contact angle in presence of boundaries. 
In order to calculate the density profile we need to use the general expression of the 
pressure tensor ()18|) and insert it in to the mechanical equilibrium condition: 



V P(x) = 0, 



(21) 



with the appropriate boundary conditions p(— oo) = p g and p(+oo) = p\ for a planar 
interface between liquid and gas. Let us suppose that the interface develops along 
the y coordinate. Under this geometry, the pressure tensor becomes anisotropic with 
a mismatch between the transverse components, P xx (y) = P zz {y) and the normal 
component P yy (y). The condition that the interface does not move, (}2*Tj) . implies 
that the normal component remains constant and equal to the value in the bulk, 
Pbuik, throughout the interface: 



yyijj^) — Pbulk 



c 2 s p + -c 2 s g b ip'- 



+ \c%^d yy ^-^\d y ^\\ 



(22) 



The density shape is now fully determined by solving (J22J) with the requirement that 
the liquid and gas phase share the same value of the bulk pressure: 



bulk 



c Ip 9 



-C 2 s Qb^ 2 {Pg 



(23) 



By making the change of variables (dpjdy)' 



z. 



one may rewrite the mechanical 



equilibrium as an ordinary differential equation for the inter-particle potential where 
only derivatives with respect to the density p appear: 



Pi 



2 ,^62,2, &bcj 
k = C s p+ —C s i) + — 



lljj 2 d 

2^7 dp 




2 /#y i 

\dp) V 



The above differential equation can be integrate explicitly to give: 



2 ,2 



T 2dP 



(24) 



(25) 



Gbctm* „ 

with ip' = dip /dp. If the two extremes of the integral are chosen inside the bulk 
phases we get: z(p) =0 for p = pi and p = p g . It is easy to realize that in order to 
be compatible with the latter boundary conditions, we must require: 



Pa 



Pi 



P, 



bulk 



C 2 S P 



c 2 Gt 







(26) 



the two densities p\ and p g as shown bv lShan &: Chen 



which fixes, together with 
( 19941 ). 

The whole profile can be obtained by inverting (|25j) or by directly numerically 
solving (|24|) starting from the inside of one bulk phase and making a small initial 
spatial perturbation on the constant density profil e. Let us notice tha t equations 
(I22l25l26j) are different from equations (21,24,25) of IShan fc ' Chenl (|l994h because of 
the different requirements imposed here in the derivation of the pressure tensor as 



2.2 Surface tension and Contact Angle 

Following Rowlinson fc Widom ( 1982| ). we define the Liquid-Gas surface tension, 



jig, as the integral along the coordinate normal to the interface of the mismatch 
between normal and transverse components of the above tensor. More precisely, 
assuming that the only dependence is in the y coordinate, the local increment of the 
surface tension is: 

Pyy-P* X = d -^ = -\ctg b \d y ^\\ (27) 

Upon integration from — oo to +00, the surface tension is readily calculated: 

1, 



li g = --4Q b J \dyWdy. (28) 

This implies the existence of transverse stresses that would result in pressure drop 
for spherical interfaces and in the most general Laplace law for a curved surface: 

** = A-k + k) (29) 



being - R1 and R?. th e local principal radii of curva ture of the surface, see lRowlinson fc Widom 



( 1982h : lDe Gennesl (l2003f ): IShan fc Che"n1 (Il994h . Let us now discuss how to general- 



ize the previous results to the case of a solid boundary In our language, mechanical 
equilibrium between multi-phase systems (say liquid, gas and solid) can be formu- 
lated as a more general problem. We want to study the mechanical equilibrium 
by imposing that the density profile must match some given value at the boundary 
position: 

Vtp (a;) = (3Q) 

ip{p{x w )) = i){p w ) 

where the solid wall is at position x w . Equation (}3T)j) is nothing but the mechanical 
equilibrium of a multiphase system in the presence of a boundary condition. The 
value of the inter-particle potential at the wall, ip(p w ), is a free parameter in the 
model, and it is not meant to be related with the "true" density of the solid phase. 
It will be used to tune different wall properties. 

Focusing, for the sake of simplicity, on two-dimensional systems (see figure IT2l in 
Appendix B) the mechanical equilibrium equation translates into 

d x P xx + d y P xy = (31) 

or equivalently {d x P xy + d y P yy = 0), whose meaning is that the flux over an arbitrary 
contour (Green's theorem) of the vector {P xx ,P xy ) is zero. With reference to figure 
IT2l in Appendix B we notice that for such a calculation we need to specify the 
pressure tensor along a solid-gas interface and liquid-solid interface. If we choose 
the rectangular contour shown in figure fT2l and impose the flux of the above vector 
exactly to zero (details are given in the Appendix B) we obtain: 



\d y ^dy- / \d y ^dy 
cos (9) = m T — (32) 



where / \d y ip\ 2 dy, I \d y tjj\ 2 dy, / \d y ip\ 2 dy indicate the positive integrals calcu- 

J si J sg Jig 

lated along the solid-liquid, solid-gas and liquid-gas interfaces. 

From the above relation ([32)1 . which fully determines the contact angle in our 
LBE scheme, one naturally extracts the definition of the surface tensions: 

To,/? = -\ c t$b I \d y 4>\ 2 dy (33) 



a,f3 



where with a, f3 we mean any two among the liquid, gas and solid phases. Notice 
that, rigorously speaking, the surface tensions between liquid-solid and between 
gas-solid are defined only modulo an additive constant, being operation ally defined 
i n ter ms of the contact angle which depends only on their difference ( De GennesI 



( 2003| )). The above definition is consistent with the requirement to have 7; s = 



when the wall has the same density of the liquid (perfect wetting). Correspondingly, 
one may also imagine a situation when the gas phase perfectly matches the wall 
properties, p g = p w , where we have ^ sg = and consequently a complete dewetting. 

2.3 Analytical and Numerical results 

From the mechanical equilibrium condition (|3()|) . we may calculate the density vari- 
ation along an interface between any two of the three phases, generalizing the calcu- 



lation made bv lShan fc Chenl f)l994r ) for the gas-liquid interface only. To accomplish 
this, we must integrate the equation for the normal component of the pressure tensor 
(j22J) imposing the boundary conditions at the solid, p(y = y w ) = p w , and p(oo) = pp 
where with (3 we denote either the liquid or the gas phase (see also section 12.41 for 
more details). In figured we show two such profiles. In our first case at the center of 
the channel there is only liquid and the gas phase cannot develop, i.e. the averaged 
density is larger than the liquid density at coexistence, while in the second case the 
averaged density is chosen such that a gas phase can develop between the liquid and 
the wall. 

Once the density profile is known , it is easy to obtain the corresponding surface 
tension by plugging the density profile into the expression (J3*3|) . For example, in 
figure 121 we show the surface tensions, •yis,1i g ,1 sg at a given temperature and 
by changing p w , i.e. by varying the wetting properties of the surface. Of course, 
only the first two will depend on p w and, accordingly to our previous discussion, 
they must satisfy 7^ = for p w = pi (perfect wetting) and ^ gs = for p w = p g . 

Given the surface tensions for a fixed temperature one may readily calculate the 
contact angle which describes the macroscopic properties of the surface as a function 
of the mesoscopic boundary condition imposed on ip(x w ) = ip(p w ). This is the first 
important methodological result of this paper, i.e. we provide a mesoscopic way 
to parametrize the hydrodynamical behavior of fluid in the proximity of a surface 
(with different contact angles) as a function of the tunable free parameter ip(p w ). 
In figure El we show the analytical results obtained by inserting the density profile 
from (|2*4"Jl into (j32H33|) . In the same figure we also show the numerical results ob- 
tained by running the LBE code and with the estimate of the contact angle obtained 
with a goniometer as shown in figure El As one can see from figure El the method is 
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Figure 1: The stationary state configuration for the density field as a function of the 
relative distance from the wall. We present the density field normalized to the center 
channel density for two different values of the average density of the system: (p) = 
1.1/3; (□) and (p) = 0.7pi (o). In both cases we integrate numerically the Lattice 
Boltzmann Equation with r = 1 in a 2d channel with two walls (bottom and lower) 
and periodic boundary condition in the stream-wise direction. The dimensions of 
the channel are L x x L y = 80 x 45 and £/& = —6.0, being Q c = —4.0 the critical point 
of the system. The liquid and gas density for this value are respectively pi = 2.65 
and p g = 0.07 in Lattice Boltzmann units. Notice that the density tends to match 
a given value at the wall (for y/L y = 0) which is different from both liquid and gas 
values, this is because of the chosen boundary condition, ip(p w = 0.5). 
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Figure 2: Surface tensions in a 3 phase system for Q b = —6.0 as a function of p w . 
The horizontal line represents the surface tension between liquid and gas (jig). The 
surface tension between solid and gas (jsg) is zero f° r Pw = Pg = 0.07 and equal to 
•fig for p w = pi = 2.65. The surface tension between liquid and solid (7^) reflects 
the opposite behavior. 
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Figure 3: Contact angle for the case Qf, = —6.0 as a function of p w . The analytical 
estimate is compared with the experimental results of Lattice Boltzmann simulations 
(+). The Lattice Boltzmann simulations are carried out in a 2d domain periodic 
along the stream- wise direction and with two walls (upper and lower wall). The 
grid mesh used is L x x L y = 80 x 45 and the relaxation time chosen is r = 1. To 
produce the steady state with a drop of liquid the system has been initialized with 
a non-homogeneous condition with a square spot of liquid on the lower wall. 



6 ~ 0°, the latter case is obtained when p w is chosen equal to the liquid density at 
the given temperature. All numerical results are obtained by using the shape fl20(l 
for the inter-particle potential. This is the first study, to our knowledge, where a 
systematic analytical procedure to derive the contact angle in Lattice Boltzmann 
models with interparticle potentials has been proposed. Other important attempts 



were already published bv lBriant et all (j2003[ ) c oncerning a la ttice transcription of 



mean field thermodynamic boundary conditions (ICahn I (jl977t )) in the framework of 



free-energy multiphase methods (see ISwift et all (J1995J)). Numerical investigations 



of th e contact angle within a LBE approach were also presented by iKwok (1200,4 



20041 ) but without any analytical control on the links with the surface tension as 



proposed here. 

The production of a rarefaction zone close to the wall is also important to de- 
termine the slippage properties in the case when a constant external pressure drop 
is applied on the system (Poiseuille flow). For instance, the formation of a gas layer 
is believed to be the most probable ca use of the apparent slip length measured in 
many experimental micro- channels fsee lTrethewav fc Meinhart ( 20021 )). The idea is 
that the rarefaction layer leads to two feedback on the bulk fluid velocity. First, it 
allows for the fluid to slide on it without touching the boundary. Second, it gives 
an effective reduction on the channel width seen by the fluid, leading to an overall 
increase of the mass throughput for a given pressure drop. Once the density profile 
is determined, one may solve the equations for the stream-wise velocity (|15j) profile 
in presence of a unitary external pressure gradient, VP ex t, arid a unitary kinematic 
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Figure 4: Left: typical density profile from which we extract the contact angle. 
The contact angle in the numeric is extracted by estimating the angle made by the 
surface where the density profile is equal to the average value (pi — p g )/2 being 
pi = 2.65 and p g = 0.07 for the present case with = —6.0. The computational 
parameters are the same used for figure H3 To produce the three phases equilibrium 
the system has been initialized to a non-homogeneous condition with a a square 
spot of liquid in contact with the lower wall. Right: a vertical snapshot of the left 
panel for x = 45. 



viscosity, u, in the Stokes approximation: 

V P 

d y {pd y u x ) = J xt = 1. (34) 

In figure 03 we show a few examples of the velocity profile at changing the contact 
angle. All data have been obtained by imposing that the flow has the fluid density 
(pi) at the center of the channel and that it ends with a density (p w ) at the boundary. 
Notice that, by increasing the hydrophobic property of the surface (increasing 8), the 
rarefaction layer becomes more and more singular, i.e. the velocity profile develops 
higher and higher shear rates at the boundary. 

To quantify this effect, one usually introduces an apparent slip length, X s , defined 
as the length were the bulk, Poiseuille-like, velocity profile extrapolates to zero away 
from the wall, see figure |H1 A simple phenomenolo gical estima t e of A . 9 for small 
density variations has been proposed previously in Benzi et al. ( 20061 ) where the 



authors show that the apparent slip length can be estimated as: 

K ~ —6 (35) 

Pw 

where with Ap we mean the density jump from the bulk to the wall va lues and with 
5 is pr oportional to the rarefaction layer, for a similar study see also Harting et al 



f|2005h . Another important quantity is the mass flow rate gain, i.e. the ratio be- 



tween the actual mass flow rate and the mass flow rate obtained assuming a perfect 
Poiseuille-like profile in the whole channel section: 




Figure 5: Velocity profiles for pressure driven Stokes flow. Profiles are obtained 
integrating numerically the Stokes equation for unit pressure drop and kinematic 
viscosity (|34|) in a homogeneous (in the stream wise direction) channel with an 
height L = 50 grid points. The profiles are then normalized with respect to the 
center channel velocity of the Poiseuille flow (u p c ) for the same geometry. From top 
to bottom we show different values of p w corresponding to different contact angles: 
p w = 0.13 (6 = 166°), p w = 0.27 (0 = 145°), p w = 0.49 (9 = 110°). 



where L is the total channel height. The mass flow rate is the only quantity which 
can be easily measured in microdevices experiments. It therefore plays an important 
role as a benchmark for many modeling methods. In figure El we show the results for 
the mass flow rate gain of our mesoscopic model as a function of the contact angle. 
Notice that one can easily gain a factor of the order of 4 4- 5 times larger for the 
case of high hydrophobic boundaries. In the same figure (right panel) we also show 
the phenomenological estimate for the mass flow rate gain obtained by using the 
expression (pTH)l with Ap/p = (pi — p w )/p w and 5 estimated as the distance (starting 
from the wall) needed to reach the 98% of the center channel density (again p{). 



2.4 Consistency with Thermodynamics 

In this section, we want to discuss a few issues on the thermodynamics consistency of 
the SC model here studied. Let us first start from the simple case of a liquid-gas in- 
terface in absence of boundaries. Up to now, in order to get the whole density profile, 
we have used the mechanical equilibrium condition supplied with the appropriate 
boundary conditions for the densities in the two phases, see eqs. (|22H26|I . Next, we 
investigate the relation between this mechanical condition and the thermodynamic 
conditions on both the density values and on the density profile (i.e. on the surface 
tension expression). The Maxwell constr uction which determines the thermody- 
namics consistency (see also ICallenl (|l985h in the (P, V) dia gram of the liquid-gas 



interfaces is built in terms of the requirement that: j^'[Pb u ik — Pb(p)]dV = where 
V oc 1 1 o. It is easv to realize that this condition leads to the following' integral 




Figure 6: Slip length for the profiles in the Stokes approximation. In the left panel 
we show an extrapolation of the profiles showed in the previous figure. The slip 
length is obtained extrapolating to zero the parabolic profile obtained as a best 
fit in the bulk region of the channel. In the right panel we show the slip length 
estimated indirectly from the mass flow rate gain with respect to the Poiseuille 
profile for different contact angles (o) and compare the with the phenomenological 
result (jH5j) . solid line. The good agreement is obtained with a prefactor in (JHHj) fixed 
to be 1.2. 



constraint for the two densities, pz,p 9 , at phase coexistence: 
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bulk 
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v> 2 (p) 



p 1 



-dp = 



(37) 



which coincides with (|26|) only if ip(p) oc p. This is the indication that the SC choice 
ip{p) = (1 — exp (— p/po)) is thermodynamically inconsistent, although one may 
notice that the discrepancy are extremely small in all realistic cases. For instance in 
figure we show the good agreement of the liquid and gas density values obtained 
by using the mechanical equilibrium equation f!26[) and the Maxwell construction 
(|3T|) for different temperatures (Gb)- 

Another problem using the SC approach to describe the interface shape consists 
in the expression for the surface tension (|28|) which is proportional to {d y ip) 2 instead 
of being proportiona l toJd y p) 2 as required by thermodynamical arguments (see 
Rowlinson fc Widom ( 1982h ). Also in this case, however, the situation is quite 
encouraging. In fact, let us notice that starting from (|20|) and the expression of the 
bulk pressure (JUJ), the critical point of the system is identified by the relations: 



dPb(p) 
dp 

These are equivalent to 

m 2 = 

Using the functional form ip(p) = 



0; 



d 2 P b {p) 
dp 2 



0. 
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(3f 



H = ~~p~ • (39) 
yb 

'1 — exp(— p/po)) it is readily checked that: 
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Figure 7: Results obtained from the mechanical stability conditions (|26|). □, and 
thermodynamic coexistence curve f)37j) . o, for liquid and gas densities. The results of 
numerical simulation are also reported (x). The numerical simulations are carried 
out in a fully periodic 2d setup with grid mesh L x x L y = 90 x 90 and a relaxation 
time t = 1. The system is then initiated to have a flat interface along the x direction. 
The top branch refer to the liquid density, p h while the bottom branch to the gas 
density, p g . 



Now, using (jHSESUl) we obtain 



m 2 = -W? = -W = -4r- (41) 
Po PoGb 

Since tp' = dip /dp, this suggests the following scaling relation 

l^| 2 ~-^H<V| 2 (42) 
which in turn would imply that the correct matching 

li 9 = —AG, / \d y iP\ 2 dy ~ -i- / \d y p\ 2 dy. (43) 

The previous argument, although exact at the critical point (Qi = Q c ), is semi- 
quant it at ively correct for \Qb\ > \Q C \. In fact, in figure |H1 we show the comparison 
between the surface tension measured in our LBE approach and the one defined by 
(J43|) . The agreement is quite satisfactory for all values of temperature, showing that 
the model is not far from being consistent also on that side. 

2.5 Two-phase mesoscopic model vs. Navier-Stokes equa- 
tions 

In this section we discuss the interplay between the mesoscopic two-phase approach 
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Figure 8: Surface tension as a function of the interaction parameter Qf,. The surface 
tension calculated using both expressions in equation the case with d y p (o) and 
the case with d y ip (□). To show the importance of the normalization factor 
we also show the surface tension calculated by simply replacing d y ip — > d y p in (J2%|) . 
(filled squares). 




description at the hydrodynamical le yel with ad-hoc bound ary conditions. A similar 
study has already been proposed by Prieziev et al. ( 2005| ) where a comparison be- 
tween microscopic MD and macroscopic Stokes flows with slip boundary conditions 
was presented. In our case, we focus on two analytical results obtained for Navier- 
Stokes equations in a channel with longitudinal or transverse (with respect to the 
flow direction) free-slip strips, i.e. the case when some inhomogeneous material is 
deposited at the surface so as to drastically change the boundary conditions of the 
Navier-Stokes field from no-slip, uu = 0, to free-slip, d n uu = 0. Here uu stands for 
the velocity component along the surface and d n is the derivative along the surface 
normal direction. 

This problem can be attacked at a purely hydrodynamical level, bypassing com- 
pletely the chemical-physical reactions at the wall which generates these two dif- 
ferent boundary conditions and focusing only on the bulk liquid incompressible 
phase. Using c onfer mal mapping on plane surfaces made up of longitudinal strips, 
( Philip! ( 1972alfbl )) proposed an exact analytical solution for the Stokes problem 
and only very re c ently a similar calculation has been proposed for transverse strips 
( Lauga fc Stolel lj200ah ). In a previous paper we have shown that the same ana- 
lytical results can be obtained within t he realm of LBE f or single flows, but with 
properly modified boundary conditions ( Benzi et al. ( 200(t| )). The motivations were 
very close to the Navier-Stokes continuum approach: neglect microscopic details 
very close to the boundary, including the possible presence of a rarefaction layer 
and try to mimic the bulk profile by renormalizing the fluid-wall interactions into 
a suitable wall function. Next we sho w th at indeed one c an recover both the ana- 
lytical results of lLauga fc Stone! l|2003h and lPhiliol (|l972allbh and the numeric of the 
single phase LBE, bv using the present two-phase model. In this way, we fill the 
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Figure 9: Schematic configurations for the simulations of laminar flows with bound- 
ary conditions of free shear and no slip. The alternating pattern of longitudinal 
(left) and transverse (right) strips of free shear is produced by varying the bound- 
ary density from complete wetting [p w = p h contact angle 0°, no slip) and perfect 
dewetting in a strip H (p w = p g , contact angle 180°, free shear). With these setups 
the 'slip' probability is £ = H/L being L = L x (L y ) for transverse (longitudinal) 
strips. All the details of numerical simulations are give in figure 



chemically active wall and we are able to follow also the dynamics inside this latter 
layer. To accomplish this goal, we performed LBE simulations of the SC model in a 
channel where ip{p w ) had a periodic structure with alternating strips of hydrophobic 
(Pw — Pg: local contact angle 180°) and hydrophilic (p w = pi, local contact angle 
0°) materials. In figure ITU1 (left panel) we show the analytical results for the effec- 
tive slip length X s obtained (for both longitudinal and transverse strips) by means 
of the NS approach and by the present LBE-SC method. The agreement between 
the two is very good. In the right panel, in order to capture the meaning of local 
boundary conditions obtained using alternating patterns of wetting and non wetting 
strips, we compare our multiphase mesoscopic approach with the integration of the 
incompress ible lattice Boltzm ann Equation with alternating strips of no-slip and 
free shear ( Benzi et al. ( 2006^ ). As we can see, while the non wetting strips simply 
reduce to a parabolic Poiseuille profile, the presence of the non wetting strips acts 
as a free shear zone where, due to the presence of the gas, the liquid can slide away 
producing the local free-shear condition. 



2.6 Interaction with the wall: the role of the density gradi- 
ent 

Up to now we have modeled the wall properties using a single parameter, ip(p w ), 
which fixes the density of the flow at the boundary. These parameters can be 
intended as the counterpart of the interaction energy between solid-liquid molecules 
in an atomistic approach. In principle, also the typical interaction length may play 
a crucial role, the existence of a characteristic distance is the result of the interplay 
between attractive and repulsive terms in the Lennard- Jones potential. In order to 
mimic this effect, it has been proposed to enrich the mesoscopic LBE by assuming 
that the interaction with the wall is sunnlemented bv another external force F... 
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Figure 10: Effective slip lenght as a function of the free shear percentage (£) in 
multiphase approach with alternating strips of wetting and non-wetting properties. 
The numerical simulations are carried out in 3D setups with L x x L y x L z = 64 x 
64 x 200 with periodic boundary conditions in the stream-wise (x) and span-wise 
directions (y). The interaction parameter is Qb = —6.0 and the relaxation time 
is t = 1. The corresponding liquid and gas densities (pi = 2.65, p g = 0.07) are 
used to mimic no slip and free shear respectively. Then, the steady state effective 
slip lengths for transverse (□) and longitudinal (o) strips normalized to the pattern 
dimension (L x for transverse and L y for longitudinal) are compared with exact 
analytical estimates given for stokes flow with alternating non slip and free shear 
(lines). Right panel: The velocity profile in the inlet of the channel and in the middle 
of a free shear strip from the integration of the incompressible lattice Boltzmann 
equation with mixed boundary conditions (straight lines) are compared with the 
results of the mesoscopic multiphase approach (□ in the inlet and o in the middle of 
free shear strip). The incompressible Lattice Boltzmann equation has been used with 
the same relaxation time and an average density equal to the one of the mesoscopic 
approach. Both simulation have been forced with a constant pressure gradient so as 
to reproduce a center channel velocity equal to u c = 0.04 in LB units. 




normal to the wall and decaying exponentially ( Kwok ( 2004 ); Sullivan ( 198 lh ): 



F w (x) = Q w p[x)e 



-\x-x w \/ri 



(44) 



where x w is a vector running along the wall location and 77 the typi cal length-scale of 
the fluid- wall interaction, also known as the Kac range parameter ( Sullivan ( 198 lh ). 
Equation (JUj) has been previously used in literature in connection with a slightly 
different LBE scheme, to show numerically how the wetting angle depends on the 
ratio G w I Qh in the presence of phase coexistence between vapor and liquid (jKwok 
(l2004h ). 

The introduction of an external force, exponentially damped in the bulk of the 
flow, allows the model to control also the gradient of the density profile at the wall. 
This is a two parameter model able to fit with high accuracy the value of the density 
at the wall, through p w , and the derivative of the density at the wall, d y p, through 
the Q w term in (|44j) . 

For the case with Q w , one must modify the structure of the pressure tensor as 
follows 



P 
1 y 



c\p + \c\Q h ^ + ~c\MM> + \g b c 4 s \ V^| 2 



'1.1 



p{s)e- sh di 



(45) 



The off diagonal term, P xy , remains the same, while the mismatch between the 
pressure tensor parallel to the interface, P xx and the pressure term perpendicular to 
the interface, P yy is changed to: 



P = P 

1 xx yy 



-^Gbipdydyip - Q u 



p(s)e- s/v ds. 



(46) 



This implies that the previous estimate 
by: 



of the contact angle must be replaced 



-^T p(s)e-^ds) dy - [ ( |<9^| 2 - ^ [* p{s)e-^ds ) dy 



GbCg Jo J J si V GbCg Jo 

cos (9) = — j 

/ \d y M 2 dy 

Jig 

(47) 

The physics at the boundary now depends on two parameters, which may change 
the contact angle and the density profile independently Unfortunately it is very 
difficult to make any quantitative measurement of the density profile in experiments 
or MD simulations. It is therefore difficult to asses in a systematic way the potential 
of the LBE model. In figure ^2 we show the contact angle that can be measured 
and calculated using (jlTjl in our LBE approach. In the same figure (right panel) we 
also show the variation in the density profile for a few values of the couple (Q w , p w ) 
leading to the same contact angle. As one can see, the model is very sensitive 
to different choices of the two free parameters, which makes it potentially useful 
to describe physical situations with large variations in the densitv profiles in the 
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Figure 11: Contact angle for a two parameter model. Left panel: for a given value 
of p w we show the contact angle for different values of Q w : p w = 1.3 (o), p w = 0.6 
(□). The analytical estimate using the mechanical stability is compared with lattice 
Boltzmann simulations. The details of the numerical simulations are the same of 
figure 01 with the introduction of a wall-fluid force as in (|44[). Left panel: for the 
same contact angle (130°) we report different profiles obtained for different choices 
of the parameters (p w ,G w ). In particular the value of p w chosen is: p w = 1.3(o) 
p w = 0.85(D) p w = 0.3(x) and the corresponding value of Q w has been chosen to 
reproduce the same contact angle of 130° estimated from (|47|). 

3 Conclusions and perspectives 

We have presented a mesoscopic model, based on the Boltzmann Equations, for the 
interaction between a solid wall and a non-ideal fluid. The model is an extension of 
the SC model for dense fluids in unbounded domains. We have first derived an ana- 
lytical expression for the contact angle and for the surface energy between any two of 
the liquid, solid and vapor phase, by introducing a parameter, 4>{p w ) which fixes the 
density value of the Boltzmann molecules at the solid wall. We have shown how in 
this way one can cover the whole range of contact angles, from a super-hydrophobic 
surface 9 ~ 180° to a condition of perfect wettability 9 ~ 0°. Concerning the thermo- 
dynamic consistency of the model we have shown that although not formally verified 
it does not introduce any systematic error in the results for the surface tension and 
for the liquid-gas density variations at changing the system temperature (Maxwell 
construction). We have discussed the connection between the rarefaction layer in 
the proximity of the wall and the production of an apparent slip phenomenon. We 
have also presented a comparison between the results obtained within the realm of 
our two-phase mesoscopic model and some analytical expression for the slip length 
calculated within a Navier-Stokes approach of a single phase fluid with alternating 
boundary conditions formed of free shear and no-slip strips. We have shown that the 
LBE approach is able to reproduce quantitatively the analytical results, supporting 
the statement that slip phenomena at the macroscopic Navier-Stokes level can be 
interpreted by an apparent slip induced by gas accumulation close to the boundary. 
We have also studied the effects on the contact angle and on the density profiles of 
the use of a second free parameter, Q w connected to the wall-force decaying expo- 



if one wants to control both the value of the density at the wall and its gradient. 
An important outcome of this study is the p ossibility to ful l y inte grate the two-phase 
approach with c omplex geometries (see also Verberg et al. ( 2004j )). For instance, re- 
cent MD studies ICottin-Bizonne et al. ( 2004j ) have demonstrated the existence of a 
wetting/dewetting transition in micro-channel with grooves. The effect is driven by 
capillarity forces which may expel the liquid out of the corrugation leading to an 
increase of the effective slip length and of the mass flow rate. The mesoscopic model 
here presented is fully capable to reproduce the same effects. Moreover, it may 
also implement more complex surface patterns (different corrugations) then those 
possible with MD simulations maintaining a control on the spatial and temporal 
fluctuation s at the macrosco pic scale. Results on this direction will be reported 



elsewhere ( Benzi et al. ( 2006| )). 
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4 Appendix A 

In this appendix we detail the calculation needed to perform the continuum limit for 
the nonideal pressure term. Let us first notice that the value of At just enters the 
discretized LBE equation (J 14)) as a normalizing factor with respect to the relaxation 
time r. Changing At just means to redefine r. 

Now, if we start from the expression for the i component of the interparticle 
force: 

Fi(x, t) = -g b ip(x) w ^( x + c a At, t)c\ (48) 

a 

assuming a stationary state {F(x, t) = F(x)) and Taylor expanding up to the third 
order we obtain: 

Fi(x) = -Gbtp(x) w a c l J)(x) + At Wacicidjtpix) + 

a a,j 
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If we choose the weights w a to satisfy the following tensor relations 

a 

a 
a 

'Yw a c % a c j a c k a c l a = Cg(SijSki + Sikdji + S u 5 jk ) 




we end up with the following expression for the z-component of F: 



At (At) 3 



(49) 



being A the Laplacian operator. The above expression for the interparticle force can 
be easily translated into an excess pressure with respect to the ideal gas expression 
(c 2 s p) using the definition: 



-djPij - d t (c 2 sP ) = Fi 
and therefore we end up with a pressure tensor of the form: 
At ^ ,„ (At) 3 g b ci(At) 3 
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(50) 



c^Qbdii/jdjifj. 



(51) 

An important remark is now in order. The continuum expression just derived de- 
pends explicitly on the interparticle potential range, here of the order of c s At. The 
limit At — > would therefore imply a vanishing interaction range, i.e. the ideal gas 



limit, Pij 



c 2 sP 5. 



5 Appendix B 

In this appendix we will detail the calculation of the contact angle used throughout 
the text. Starting from the pressure tensor (for At = 1) 



% - ^Gbdi^djijj. (52) 



in the coordinate system of figure (|T2"j) we make the following change of variables 

f t. 1 = r\t\( 0) -4- ii ms( f)) 



being y' aligned with the separation line between liquid and gas. 

d_ _ dx^_d_ _|_ cV _d_ _ s [ n (Q)J> cos(8) 9 
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If we assume only the x' dependence we can write 



(54) 



ip(x,y) = ip(x') (55) 



£ = sin (0)^7 

8 #__ A (56) 



which clearly imply 



and 
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sin (0) 

= sin (e)dx' (58) 



dx — > — — t-t- / dx' (57) 



(d x tl)){dyijj)dx = , / dx' sin(0) cos (0)(d x >ip)(d x ^) = / dx' cos (0)(d x ^)(d x 'ip) 
sm (0) J y 

(59) 

By imposing a zero flux out of the boundaries of the vector (P^, P X3/ ) we obtain 
the mechanical definition of the contact angle. To this purpose, let us notice that 
based on the definition of the pressure tensor, we can write: 

c 4 

Pxx = P yy + -fGbdy^dy^fj (60) 

Pxy = jGb^dxdyij (61) 

with Pyy constant along the solid-liquid and solid-gas interface (segments A, B 
of figure H2J). This, together with (|53|) implies 
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For the case with Q w , the previously used pressure tensor must be slightly modified 
due to the presence of a normal (say along the y direction) force between the fluid 
and the wall: 



1 
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<p + A^r + A&^ + -g b ct\ v# 



i r v 

^ij-AQbdi^djij+Siyg^j / p( 



s)e 



Then relations (|60|) and (|61|) should be slightly modified since P xy gives the same 
result but for P xx we have 

P xx = P yy ~ $g b d y i>d y ^ - G w I" p(s)e- s/ri ds (65) 
^ Jo 

that imply the following estimate for the contact angle: 




/ \d x ^\ 2 dx' 

Jig 

(66) 

Notice that the importance of wall effects appears only in relation to bulk terms in 
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